Effective s- and p-Wave Contact Interactions 
in Trapped Degenerate Fermi Gases 
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The structure and stabihty of dilute degenerate Fermi gases trapped in an external potential 
is discussed with special emphasis on the influence of s- and p-wave interactions. In a first 
step an Effective Contact Interaction for all partial waves is derived, which reproduces the 
energy spectrum of the full potential within a mean-field model space. Using the s- and p- 
wave part the energy density of the multi-component Fermi gas is calculated in Thomas-Fermi 
approximation. On this basis the stability of the one- and two-component Fermi gas against 
mean-field induced collapse is investigated. Explicit stability conditions in terms of density and 
total particle number are given. For the single-component system attractive p-wave interactions 
limit the density of the gas. In the two-component case a subtle competition of s- and p-wave 
interactions occurs and gives rise to a rich variety of phenomena. A repulsive p-wave part, 
for example, can stabilize a two-component system that would otherwise collapse due to an 
attractive s-wave interaction. It is concluded that the p-wave interaction may have important 
influence on the structure of degenerate Fermi gases and should not be discarded from the 
outset. 

PACS numbers: 03.75.Fi, 32.80.Pj, 34.20.Cf 



I. INTRODUCTION 

The achievement of Bose-Einstein condensation in 
trapped dilute gases of bosonic atoms triggered a 
wide spread interest in the field of ultracold atomic gases. 
Meanwhile these systems appear as a unique lab for the 
study of all kinds of fundamental quantum phenomena. 

After a series of excellent experiments on bosonic sys- 
tems the question arises, whether a dilute gas of fermionic 
atoms can also be cooled to temperatures where quan- 
tum effects dominate. In 1999 the group of Deborah 
Jin managed to cool a sample of typically 10^ fermionic 
■*'^K atoms to temperatures significantly below the Fermi 
energy of the system Q; in recent experiments they 
achieved a temperature corresponding to one fifth of the 
Fermi energy Q . In this temperature regime the system 
can be described as a degenerate Fermi gas, where the 
majority of the atoms successively fills the lowest avail- 
able one-body states according to the Pauh principle. 

One of the goals of the investigations on dilute and ul- 
tracold Fermi gases is the observation of Cooper pairing 
and the transition to a superfluid state. The transition 
temperature depends on the density and on the strength 
of the attractive interaction that is necessary for the for- 
mation of Cooper pairs |5|. In order to increase the 
transition temperature one may increase the density or 
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the interaction strength, where the latter seems to be 
more promising. An atomic species favored for the exper- 
imental observation of a BCS transition is ^Li due to its 
large s-wave scattering length of oq « — 2160aB But 
also ''"K, which shows a rather small natural scattering 
length, is a possible candidate for Cooper pair formation 
Q , because a simultaneous s- and p-wave Feshbach reso- 
nance is predicted which allows tuning of the s- and 
p-wave scattering lengths over a very wide range. 

A serious constraint on the way towards a superfiuid 
Fermi gas is the mechanical stability of the progenitor, 
i.e., the normal degenerate Fermi gas. As was experi- 
mentally demonstrated for the bosonic *^Rb system [||, 
an attractive interaction between the atoms leads to a 
mean-field instability of the trapped gas if the density 
exceeds a critical value. A similar collapse is expected in 
fermionic systems with attractive interactions. In con- 
trast to the bosonic systems p-wave interactions con- 
tribute in a Fermi gas and may have strong influence 
on the stability of the system |^ . 

In the following we address the question of the stabil- 
ity of degenerate one- and two-component Fermi gases 
in presence of s- and p-wave interactions within a simple 
and transparent model. In Section || we derive an Effec- 
tive Contact Interaction (ECI) for all partial waves that 
is suited for a mean-field treatment of the many-body 
problem. Using the s- and p-w ave part of this interac- 
tion we construct in Section III the energy functional of 
a trapped multi-component Fermi gas in Thomas-Fermi 
approximation. From that the ground state density pro- 
file can be determined by functional variation. In Sec- 
tion IV we discuss the structure and stability of single- 
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component Fermi gases, where only the p-wave interac- 
tion contributes according to the Pauh principle. Sec- 
tion ^ deals with two-component Fermi gases where s- 
and p-wave interactions are present and lead to a subtle 
dependence of the stability on the two scattering lengths. 



II. EFFECTIVE CONTACT INTERACTION 

A. Concept 

The approximate solution of the many-body prob- 
lem in a restricted low-momentum sub-space of the full 
Hilbert space faces a fundamental problem: Many realis- 
tic two-body interactions — like van-der-Waals-type in- 
teractions between atoms or the interactions between two 
nucleons in an atomic nucleus — exhibit a strong short- 
range repulsion. This repulsion generates particular 
short-range correlations in the many-body state, which 
cannot be described within a low-momentum model- 
space . Therefore short-range correlations inhibit the 
use of a realistic two-body interactions in the framework 
of a naive mean-field model. One way to overcome this 
problem is to replace the original interaction by a suit- 
able effective interaction. How this effective interaction 
has to be constructed from the original potential depends 
on the properties of the actual physical system under in- 
vestigation. 

Cold and very dilute quantum gases allow an effective 
interaction of simple structure. The typical wave length 
of the relative motion of two particles is always large com- 
pared to the range of the interaction. Therefore the par- 
ticles experience only an "averaged" two-body potential 
and do not probe the detailed radial dependence. More- 
over the gases are in a metastable not self-bound state 
that is kept together by the external trapping potential. 
Thus the bound states of the two-body potential are not 
populated and have only indirect influence. We make 
use of these facts and replace the original potential by a 
contact interaction. The strength of the contact terms is 
related to the properties of the original potential. 

We establish this relation via the energy spectrum of 
the two-body system. Consider a two-body problem with 
an auxiliary boundary condition at some large radius 
such that the spectrum is discrete even for positive ener- 
gies. For noninteracting particles this spectrum shows a 
sequence of levels at positive energies Eni which are la- 
beled by an angular momentum quantum number I and 
a radial quantum number n. A schematic sketch of this 
spectrum for some fixed I is shown in Figure |^. If we 
switch on an interaction between the particles two things 
happen: A number nf of bound two-body states with rel- 
ative angular momentum I may appear at negative en- 
ergies and the levels at positive energies Eni are shifted 
compared to the noninteracting spectrum. For our pur- 
pose this energy-shift AEni — Eni ~ Eni is defined be- 
tween the positive energy states only. The lowest level 
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FIG. 1: Schematic comparison of the free two-body energy 
spectrum Eni with the energy spectrum Eni in presence of 
a two-body interaction for given relative angular momentum 
I. The interaction generates nj' bound states and shifts the 
positive energy states by AU^i. 



of positive energy is labeled with the quantum number 
n = 0. 

The energy-shift is the relevant "averaged" property 
of the two-body potential that the effective interaction 
should reproduce. We require that the expectation value 
of the effective interaction calculated with eigenstates 
\nlrri) of the noninteracting two-body problem equal the 
energy shift A£^„i induced by the original potential 



{nlm\ v""^ \nlm) = Eni - E^i = A£;„; 



(1) 



Due to the angular momentum dependence of the energy- 
shifts we have to construct the effective contact interac- 
tion as a sum of contact interaction for each partial wave. 
The strength of each contact term is fixed by this condi- 
tion. 



B. Construction of the ECI 

The construction of the Effective Contact Interaction 
(ECI) is organized in two steps: Firstly, the two-body 
energy-shift induced by the original potential is evalu- 
ated. Secondly, the operator of the Effective Contact 
Interaction is formulated for each partial wave and the 
interaction strengths are fixed by condition (^. 



1 Energy Shift 

Assume a system of two particles with reduced mass 
/i = TOiTO2/(™i + ^2)- Their wave function can be de- 
composed into a center of mass and a relative wave func- 
tion, where the latter separates into a radial and a angu- 
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lar component because of rotational symmetry 

{r\nlm) ^ Rni{r)Yi,n{n). (2) 

The radial wave functions Rniir) of the noninteracting 
two-body system are solutions of the Schrodinger equa- 
tion 



1 1(1 + 1) ^ ^ 



Rniir)^0 (3) 



(we use units with h ^ I). We require the radial wave 
function to vanish at some arbitrary but large radius A 



Rnl{A) = 0. 



(4) 



This auxiliary boundary condition leads to a discrete en- 
ergy spectrum, which is needed to enumerate the energy 
levels and to evaluate the energy shift. 

In the noninteracting case the solution of the radial 
Schrodinger equation (^) is given by a spherical Bessel 
function ji (x) 



Rni{r) = Anijiiqnir), 



(5) 



where q„i denotes the relative momentum of the two par- 
ticles and Ani a normalization constant. The discrete 
momenta qni are determined by the boundary condition 
(^) and thus are related the to the zeros of the Bessel 
function. Since the radius A can be chosen arbitrary 
large it is sufficient to use the asymptotic expansions of 
the spherical Bessel and Neumann functions 



ni (x) 
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7r//2), 
- 7r//2). 



(6) 



Evaluating the boundary condition (Q) with the asymp- 
totic form of the Bessel function we obtain for the possi- 
ble relative momenta 



g„(A = 7r(n + Z/2). 



(7) 



Accordingly the two-body energy spectrum in the nonin- 
teracting case is given by 



(8) 



The normalization constant can be determined explicitly 



a: 



V A^ 
dr r^jfiqnir) = — jf+i(g„,A). (9) 



Inserting the asymptotic form of the Bessel function we 
finally get 



qniA:»i A 



(10) 



In presence of a two-body potential v{r) of finite range 
A the solution Rni (r) of Schrodinger equation 



1 Ijl + l) 

-r H h z/i[w(r 



r dr^ 



Enl 



Rniir) = 

(11) 



outside the range of the potential (r > A) is given by a 
superposition of spherical Bessel and Neumann functions 

R„l{r) = Ani [jiiqnir) - tan77/(g„;) «;(?„; r)] , (12) 

where r]i{q) denotes the phase shift of the potential for 
the ^-th partial wave. The bar distinguishes quantities 
in presence of the interaction from those in the nonin- 
teracting case. Imposing the boundary condition (|^) we 
get the following implicit equation for the momenta qni 
in the interacting case 



ni{qniA) 



tan r]i{qr, 



(13) 



Expressing the Bessel and Neumann functions by their 
asymptotic expansion (||) this reduces to 



tan(g„;A - 7rZ/2) = tanr/;(g„;) 



(14) 



In order to associate the energy levels with same quantum 
number n in the interacting and noninteracting case as 
described above, the lowest positive-energy state should 
be labeled with the quantum number n — 0. To achieve 
this for a potential with n]' bound two-body states with 
angular momentum quantum number / we add the phase 
7r(n -|- n^) to the argument on the r.h.s., i.e., the bound 
states contribute according to Levinson's theorem. The 
momenta g„; in presence of the interaction are thus de- 
termined by the equation 



qni 



A: 



-m{qni) + TT{n + n]; + 1/2). 



(15) 



The momentum shift Ag„; induced by the interaction 
is obtained by subtracting the momenta of the noninter- 
acting system (^ from equation ( [l5| ) 



Ag„iA = {qni - qni)^ 



nl) 



(16) 



'nl) 



Here fji {q) — i]i (g) — nn]' denotes the phase shift reduced 
by the contribution of the bound states. In a final step 
we expand the phase shifts around the momenta qni of 
the noninteracting system 



miqnl) = miqnl) + T)'i{qnl)^qnl + 



(17) 



Already the term linear in Aq„; can be neglected in good 
approximation because the momentum shift is of the or- 
der 1/A according to equation (|l^). Thus we retain the 
following simple expression for the momentum shift 



Aq„iA = -fiilqn 



(18) 
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The shift of the energy levels of the interacting two- 
body system compared to the noninteracting spectrum 
is connected with the momentum shift by 



Enl 



qni 



q-ai 



(19) 



The term quadratic in the momentum shift can be ne- 
glected because Aqni/qni is small. The final expression 
for the energy shift reads 



AE„ 

Enl 



_2 fii{qni) 
A qni 



(20) 



The proportionality between energy shift and phase shifts 
is well known and was used in different applications 
before. 



2 ECI in Scattering Length Approximation 

In a second step we construct an operator form of the 
Effective Contact Interaction that obeys condition (|l]). 
For the application to ultracold dilute quantum gases the 
general form (^ of the energy shift AEni can be simpli- 
fied considerably. The relative momenta in these systems 
are extremely low, i.e., typical relative wave lengths are 
large compared to the range of the interaction. This al- 
lows an expansion of the phase shifts f}i{q) in a power 
series in relative momentum. In lowest order approxima- 
tion [qai <C 1) the phase shifts of the Z-th partial wave 
can be expressed in terms of the corresponding scattering 
length ai [20] 



m{q) _ tan77,(q) 



„2l+l 



{21 + 1) 
[(2/ + l)!!]2 



(21) 



The energy shift in this scattering length approximation 
is given by 



AE„i _ 2 {21 + 1) 
Enl A[(2/ + l)!!]2 



qni 



(22) 



Based on this form we can construct a simple operator for 
the ECI. For applications where the approximation ( |2l| ) 
is not sufficient higher order terms of the power series 
can be included successively. We will come back to this 
point in the following section. 

According to the dependence of the energy shifts on 
the angular momentum quantum number I the operator 
of the effective interaction is formulated as a sum of 



independent operators v"^^ for each partial wave 



V = 



(23) 



where 11; denotes the projection operator on the sub- 
space spanned by states of relative angular momentum 
I. We want to use a contact interaction for each par- 
tial wave. This requires nonlocal interaction operators 



beyond I = 0, i.e., derivative couplings. The simplest 
ansatz for the operator of the effective contact interac- 
tion for the l-th partial wave is 



..off 



{q-rY giS^^Hr) {f-q)' 
Jd'r\r)-^9iS^^Hrl^{r 



(24) 



Here q — ^{pi — P2) denotes the operator of the rela- 
tive momentum of two particles, f = r/r the unit vec- 
tor along the relative coordinate and 5^^^(r) = j^5{r) 
the radial component of the 3-dimensional delta function. 
The arrows above the derivatives indicate to which side 
they act. 

The interaction strength gi is a constant that contains 
the relevant information on the original two-body poten- 
tial. The connection is provided by condition (|^) via the 
shift of the energy levels with respect to the free spec- 
trum. To evaluate (^ we calculate the expectation value 
of the contact interaction (^) for the l-th partial wave 
taken with the noninteracting two-body states \nlnij 



{nl 



m V, 
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47r 

9i_ 

47r 



d'rd^'\r)\^Rni{r)YUnf 

Ql 2 



(25) 



(2^ + 1)!! 



r=0 
2 



^nl qni ' 



where we used the expansion of the radial wave-function 
around the origin 



Rni {r) = Anl 



(2^ + 1)!! 



2(2/ + 3) 



(26) 



Equating the expectation value ( pSD with the energy 
shift results in the following equation for the interaction 
strengths 
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(27) 



Inserting the normalization constant ( |1C|) and the energy 
shift in scattering length formulation (E3) gives the final 
expression for the interaction strengths 



4n {21 + 11 2/+1 



2^l (n)2 



(28) 



Together with equation ( |24| ) this defines the scattering 
length formulation of the Effective Contact Interaction. 
Notice that these expressions do not depend on the aux- 
iliary boundary condition (^, which was introduced to 
obtain a discrete energy spectrum. 

For I — 1 sometimes a gradient operator with respect 
to the relative coordinate is used instead of the radial 
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derivative (^4|). This alternative ansatz reads 
= /dV|f)V giS^^^r) V(f| . 



(29) 



It should be noted that the p-wave interaction strength 
gi in the gradient formulation is related to the interaction 
strength (|2|) by 



51 = 5i/3- 



3 Beyond Scattering Length Approximation 



(30) 



The operator of the Effective Contact Interaction 
can be generalized systematically beyond the scattering 
length approximation shown in the preceding section. A 
formal scheme emerges from the expansion of the phase 
shifts in a power series in q 



oc 



1 



iy=0 



(31) 



The momentum independent lowest order term of this 
expansion matches the scattering length approximation. 
From equation ( |2l] ) the connection between the coeffi- 
cient cf^^ and the the scattering length ai becomes obvi- 
ous 



„(o) 



The coefficient c[^'' of the quadratic term of the expansion 
is connected to the so called effective range or effective 
volume of the potential. We will discuss this contribution 
in more detail later. By inserting the expansion (^ij) into 
the general formula (EO) for the energy shifts we obtain 



(2/ + 1) 



,2/+l 



(32) 



Enl 



2 

"a 



OO 

u=0 



(33) 



Following the basic concept of the ECI these energy 
shifts have to be generated by the operator of the ECI 
according to condition (^. To include the momentum 
dependence of the energy shifts we have to generalize 
the ansatz for the ECI operator compared to the sim- 
ple momentum-independent scattering-length formula- 
tion (24). The further calculation will show that 



u=0 



-f (q-f)'+2'' ,5(3)(f) (f.q)'] 



(34) 



is a proper ansatz for the effective interaction operator for 
the l-th partial wave. Besides the more complex nonlocal 
structure a set of interaction strengths gi'^^ {v = 0,1, . . .) 
for each partial wave is included. They are related to 



the coefficients cl"^ and thus correspond to the different 
powers of the momentum in equation (|3^). To employ 
condition (|l]) we calculate the expectation value of vf^ 
with nonintcracting two-body states \nlm^ 



47r 



dr'- 



Rni(r) 



r=0 



Rniir) 



r=0 



Att 



47r i^"" {21 + 1)!! 2''v\{2l + 2;^ -f 1) 



(35) 



where the full expansion of the nonintcracting radial wave 
function (p|) around r — was used |p2[ 



Rni (r) = Ani 



^^2^^1(2? + 2^ + 1)!! 



{q-nir 



(36) 



By inserting the expansion of the energy shifts ( p3| ) and 
the expectation value into condition (|l|) and comparing 
the coefficients for the different powers of the momentum 
Qni we obtain 
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(-1) 



^ +i'iTr 2'^(2;-H)!!(2? + 2;^ + l) 



2ii 



l\{l + 2v)\ 



Thus the interaction strengths g\'^^ of the general oper- 
ator form of the ECI (|3^) are proportional to the coeffi- 
cients ci^'^ of the expansion of the phase shifts (|3l|). The 
equations (Q) and ( ^ ) define the most general form of 
the Effective Contact Interaction. 

For the application to dilute degenerate Fermi gases we 
will use the ECI up to quadratic terms in the momentum, 
i.e., we include the scattering length term of the s- and p- 
wave part as well as the s-wave effective range correction. 
At this point we have to discuss the connection between 
the quadratic term of the expansion (|3^) and the usual 
effective range theory. For the s-wave phase-shifts the 
effective range expansion reads 



(37) 



q cot 770(g) 



1 1 

ao 



(38) 



where is the effective range of the potential. If we 
convert this into an expression for 770 (g) / q and expand in 
q we obtain 



-ao - 60 q -I- 



(39) 



with an effective volume &o that depends on the scatter- 
ing length flo and the effective range ro 



(40) 



Rather than using this relation we will adjust 60 in order 
to get the best representation of the phase shift fio{q)/q 
with the truncated expansion (39). 
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Finally, inserting Cg^^ = —60 into equation ( |37| ) gives 
an expression for the interaction strength of the s-wave 
effective range term 



(1) 12^. 



(41) 



C. Example: Square- Well Potential 

To illustrate the concept of the Effective Contact In- 
teraction we use the simple toy-problem of two particles 
interacting by an attractive square-well potential of ra- 
dius A and depth — Vq. 

First we look at typical wave functions to which the 
idea of the ECI applies. The major condition is that the 
typical wavelength of the relative motion is large com- 
pared to the range of the interaction. This is ensured 
by choosing the radius A of the square- well much smaller 
than the radius A associated with the boundary condition 
(^); in the following we use A = 0.01 A. Figure ^ shows 
the radial wave function of the lowest I — state with 
positive energy for different potential depths Vq. Out- 
side the potential the structure of the wave functions is 
very similar for the different interaction strengths. Only 
the wavelength is changed slightly due to the different 
matching to the wave function in the interior (see in- 
sert of Figure This change of the relative momentum 
translates immediately into an energy shift. From that 
picture the connection between energy shifts and phase 
shifts is evident. 

A second point becomes clear from this simple exam- 
ple: The detailed structure of the radial dependence of 
the potential or the number of bound states is irrelevant 
for the energy shift, only the phase shifts 'fii{q) matter. 
The insert in Figure ^ shows that the wave functions be- 
have very different within the range A due to the different 
potential depths. Moreover the potentials have a differ- 
ent number of bound states, e.g., the thick solid curve is 
associated with a potential with one bound state but zero 
phase shift. Hence the behavior outside the potential is 
identical to the noninteracting case (thin solid curve) and 
the energy shift is zero. 

Next we investigate the dependence of the / = en- 
ergy shifts on the strength of the attractive potential. 
Figure plshows the relative energy shift AE°'^'^'^^/E ver- 
sus A\/Vo, where the radius A — O.OlA of the square-well 
potential is kept fixed and the depth Vq is increased. A 
characteristic pattern appears: In the vicinity of interac- 
tion strengths where the potential gains another bound 
state the relative energy shift assumes large positive and 
negative values. Large positive energy shifts occur for 
potentials that have a very weakly bound state, nega- 
tive energy shifts for those that have an almost bound 
state. In between these interaction strengths extended 
plateaus of nearly constant energy shift appear. Within 
the plateaus the energy shift is independent of the radial 
quantum number n of the level or the relative momen- 
tum. This is a special property of the s-wave channel; 



for higher partial waves the relative energy shift ( p2| ) is 
proportional to 5^' . Only at the edges of the plateaus a 
slight dependence on the relative momentum shows up 
(see Figure |3|). 

This structure is closely related to the behavior of the 
s-wave scattering length. The I — energy shift induced 
by the ECI in scattering length approximation ( p^ ) is 
proportional to the scattering length ao- For the square- 
well potential we get 

^^2^ with ^^i_tan(AV5)^ 
Eno A A ^ ' 

This ECI energy shift is right on top of the solid curve 
in Figure^, i.e., it agrees very well with the exact energy 
shift for the lowest positive energy state. Even for higher 
momenta the agreement is very good provided that the 
magnitude of the scattering length is not too large. Sig- 
nificant deviations occur only if momentum and scatter- 
ing length are large. 

To obtain a quantitative measure of for the applica- 
bility of the ECI in scattering length approximation we 
investigate the relative deviation [E — E°^'^^^)/E of the 
ECI energy levels compared to the exact ones for the 
square-well potential. We expect that the agreement gets 
worse if either the relative momentum or if the scatter- 
ing length is large. Therefore, Figure^ shows the relative 
energy deviation for I = and I ~ 1 states as function 
of the product of momentum and scattering length, qai, 
which was assumed to be small in order to introduce the 
scattering length (pl|). The different curves correspond 
to different values of the momentum qni and were ob- 
tained by varying the depth of the square-well potential 
and thus the scattering length. 

As expected the deviation increases with increasing 
value of qai. Nevertheless, the deviation of the energy 
calculated with the ECI in scattering length approxima- 
tion is below 1% up to rather large values of qai < 1. 
If we tolerate a maximum deviation of 5% then the 
scattering length formulation can be used up to values 
qai ^ 1-5. 

It should be noted that the relative deviation of the 
general form (|2^) of the ECI in the parameter range dis- 
cussed above is below 10^*. Thus all approximations 
made to obtain equation (^0|) are valid on a high level of 
accuracy. The restrictions on the validity of the scatter- 
ing length formulation (22) originate for the replacement 
of the phase shifts by the scattering length alone, which 
is no inherent part of th ECI concept. If the simple scat- 
tering length formulation is not sufficient for a special 
application one can go beyond that. 

For example, the inclusion of effective volume correc- 
tions (see eq. (^)) improves the agreement with the 
exact energy shifts. In this way we can reduce the maxi- 
mum deviation to only 1% up to qao < 1.5. 
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rRoo{r) 
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FIG. 2: Radial wave functions rRoo{r) of the lowest I — positive energy state for different potential depths Vo of the square- 
well potential. The interaction strengths are A\/W = (thin solid), 4.49 (solid), 4.85 (dashed), and 9.5 (dotted). The radius 
of the well is A = O.OIA and marked by the gray area. The insert shows a magnification of the region around the origin. 




E-E 




1 2 

qai 



FIG. 3: Relative energy-shift AE^^o"""^ / E„o obtained from the 
exact I — solutions plotted versus the strength Xy/Vo of the 
square-well potential. The solid line gives the energy shift 
for the lowest {n = 1), the dashed line for the 10th, and the 
dotted line for 20th state of the positive energy spectrum. 
The dots mark the interaction strengths used in Figure ^ 



FIG. 4: Relative deviation of the two-body energy calculated 
with the ECI in scattering length approximation from the 
exact energy as function of qai for I = (left) and I = 1 
states (right). The curves were obtained for interactions with 
three bound states by varying the strength \\/Vo and looking 
at the energy shifts for different relative momenta g„;A ~ 20 
(solid), 40 (dashed), and 80 (dotted). 



D. ECI versus Pseudopotential 

The idea to simulate the effect of a complicated finite- 
range two-body potential by a simple s-wave contact in- 
teraction dates back to E. Fermi ||l3| and was used by sev- 
eral authors in various physical contexts. K. Huang 
and C.N. Yang |l§ generalized this idea and con- 
structed the so called pseudopotential that acts in all 
partial waves. 

The aim of the pseudopotential is (a) to generate the 
phase shifts of the original potential by a boundary con- 
dition at r = and (b) to reformulate this by an addi- 
tional inhomogenous term in the Schrodingcr equation of 
the two-body scattering problem. This additional term is 
interpreted as the pseudopotential, which can be phrased 



in the following operator form m 



pscudo 



with an interaction strength 



)2/+l 



pscudo _ 47r (/ + 1) 2/+1 

" (2Z-)-l)! 



(43) 



(44) 



for the Z-th partial wave. For this discussion we restrict 
ourselves to the scattering length approximation of the 
phase shifts. Due to the fact that the radial derivative 
acts only to the right hand side the operator of the pseu- 
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dopotential (^3|) is not hermitian. This is in contradiction 
to the basic concept of effective interactions. 

A more severe weakness shows up when we evaluate 
the energy shifts i nduc ed by the pseudopotentiaL As dis- 
cussed in section II A the expectation value of a proper 
effective interaction with eigenstates of the nonin- 

teracting two-body system should be equal to the energy 
shift induced by the original potential. For the pseudopo- 
tential (MS) we obtain an energy shift 



^pseudo 



'nl 



1 

Enl 
2 



a + 



pscudo I 7 



(45) 



A [(2? + l)!!]2 



„2i 2;+i 



This has to be compared with the full energy shift in 
scattering length approximation ( |2^ ) which is by con- 
struction reproduced by the ECI. Obviously the energy 
shift induced by the pseudopotential for states with I > 
is by a factor {l + l)/{2l + l) smaller the the energy shift of 
the original potential. Thus the pseudopotential under- 
estimates the effect of the two-body interactions beyond 
s-wave when used in a mean-field framework. For the 
widely used s-wave part the energy shifts of the pseu- 
dopotential agree with the energy shifts of the original 
potential. 

We conclude that the nonhermitian pseudopotential is 
not a proper effective interaction for a mean-field descrip- 
tion of a dilute quantum gases that goes beyond s-wave 
interactions. 



III. ENERGY FUNCTIONAL OF A TRAPPED 
MULTI-COMPONENT FERMI GAS 

A. Fundamentals 

In the following we investigate the ground state prop- 
erties of a dilute Fermi gas composed of S distinguishable 
components that is trapped in an external potential U {x) 
at temperature T = K. In present experiments |^ one 
or two components are used, which belong to the same 
atomic species but are distinguished by different projec- 
tions Mp of the total angular momentum F onto the 
direction of an external magnetic field. We distinguish 
the different components by a formal quantum number 
^ = 1, . . . , S. For simplicity we use the same mass m of 
the atoms for all components. 

We treat the many-body problem in the framework 
of density functional theory and construct an energy 
functional of the inhomogenous multi-component Fermi 
gas within a proper approximation. The ground state 
density-distribution of the many-body system is then de- 
termined by functional minimization of the energy. 

The large particle numbers of the order TV ~ 10^ allow 
the rather simple Thomas-Fermi approximation for the 
energy functional. It is assumed that the energy density 
of the inhomogenous system is described locally by the 



energy density of the corresponding homogenous system, 
and higher-order terms which include gradients of the 
density are small. To check the quality of the Thomas- 
Fermi approximation we calculated the next order gradi- 
ent corrections for a trapped noninteracting Fermi gas. 
For N — 100 particles the relative contribution of the 
gradient correction to the total energy is of the order of 
10~2; for typical particle numbers of = 10^ it drops to 
10-5 0.^ 

As starting point for the Thomas-Fermi approximation 
we calculate the energy density of the homogenous inter- 
acting multi-component Fermi gas in mean-field approx- 
imation. The basic restriction of the mean-field picture 
is that two- and many-body correlations induced by the 
interaction are not contained in the many-body state. 
Nevertheless they can be implemented implicitly by us- 
ing a proper effective interaction that is tailored for the 
model-space available. In the previous section we con- 
structed the Effective Contact Interaction especially for 
the mean-field description of dilute not self-bound quan- 
tum gases. 

A central topic of the following studies is the role of 
the interaction on the structure and stability of trapped 
degenerate Fermi gases. Our special interest concerns the 
p-wave part of the interaction, which contributes even at 
T = OK — in contrast to bosonic systems. It will turn out 
that the p-wave terms can be of substantial importance 
for the ground state properties of fermionic systems and 
should not be neglected from the outset. 

We write the Hamilton operator of the system as a sum 
of the external trapping potential U and an internal part 
Hint 



H=U + Hn 



(46) 



The internal part contains the kinetic energy and the Ef- 
fective Contact Interaction as discussed in section ||. We 
include the s-wave and p-wave terms of the ECI in scat- 
tering length formulation as well as the s-wave effective 
range correction. With equations (p^, (p8|), (|3^), and 
( ^ ) the internal Hamiltonian reads 



1 
2m 
An 
m 

m 



E 



<5(^)(r. 



(47) 



^0 E J['5^^H^»j)(nr9.j)'+h.a.] 



The summations over the particle indices i and j range 
from 1 to the total number of particles. The properties of 
the two-body interaction are parameterized by the s- and 
p-wave scattering lengths ag and oi, respectively, and by 
the s-wave effective volume 6o- In general the interaction 
parameters depend on the component quantum numbers 
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^ of the interacting particles. In order to discuss the 
basic phenomena we restrict ourselves to equal interac- 
tion parameters for all components. The generalization 
to scattering length matrices that account for the depen- 
dence on the component indices of the two interacting 
particles is straightforward. 

Experimentally each component may experience a dif- 
ferent trapping potential U(^{x). For magnetic traps this 
is due to the different magnetic momenta of the com- 
ponents, which leads to a relative shift of the trapping 
potentials for the components. Thus the operator of the 
external potential has the following form 



(48) 



where 11^ is a projection operator onto states with the 
component quantum number ^. 

Many of the results shown in the next sections do not 
depend on the actual shape of the trapping potential. If 
the shape enters explicitly we assume a deformed har- 
monic oscillator potential 



— ^ /'\2 2 I a2 2 I •v2 2\ 



(49) 



where lu — ^luiuj2Lli^ is the mean oscillator frequency 
and i — (mw)^^/^ the corresponding mean oscillator 
length, i.e., the mean width of the Gaussian single- 
particle ground state of the harmonic oscillator poten- 
tial. The deformation is parameterized by the ratios 
Xi — uji/uj, which fulfill the condition A1A2A3 — 1. 



B. Energy Density in Thomas- Fermi 
Approximation 

The calculation of the energy density functional of the 
inhomogenous interacting Fermi gas is performed in two 
steps: First we calculate the energy density of the corre- 
sponding homogenous system in mean-field approxima- 
tion. In the second step this is translated into an en- 
ergy density of the inhomogenous system by means of 
the Thomas-Fermi approximation. 

The ground state of a many-fermion system in mean- 
field approximation is given by an antisymmetrized prod- 
uct of one-body states In the case of a homogenous 
system the one-body states are eigenstate of the momen- 
tum operator with eigenvalues ki. In addition they are 
characterized by the component quantum number ^ 



(50) 



Assuming a box of volume V with periodic boundary 
conditions the spatial part of the one-body states is given 
by 



(f|fcj) 



1 



exp (i ki ■ x) 



(51) 



The energy density of the homogenous system is given 
by the expectation value of the internal part of the Hamil- 
ton operator (E^ 



^hom — {Hint)/V. 



(52) 



The calculation of the expectation values of the several 
parts of the Hamiltonian is straightforward ||l^. As a 
function of the Fermi momenta of the different com- 
ponents the energy density reads 



^hom ('^l 7 ■ ■ ■ I '^H 



+ 



5 

207r^TO ^ '^^ 



(53) 



+ 



307r3 
al + bp 



The summations run over all components ^ = 1, ... ,5. 
To avoid fractional exponents we use Fermi momenta 
rather than densities = k^/{6tt'^). 

The basic assumption of the Thomas-Fermi (or local- 
density) approximation is that the energy density of the 
inhomogenous Fermi gas is locally given by the energy 
density of the corresponding homogenous system. Thus 
the energy density of the inhomogenous system is con- 
structed from ( [53[ ) by replacing with local Fermi mo- 
menta (x) . In addition the contribution of the external 
trapping potential has to be included. This results in the 
following expression for the energy density of the trapped 
interacting multi-component Fermi gas 



^ E 



+ 



20n^m 



3 

]_ 

SOTT'^m 



E ^|(x)4(x) 



(54) 



al + bo 



E [4(f)4(x) + K|(f)4(f)]. 



The local Fermi momentum is related to the density of 
particles of the component ^ by 



1 



(55) 
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Accordingly the number of particles of component ^ is 
given by 



1 

6^ 



d^x 



(56) 



As discussed in section [I C we can reproduce the two- 
body energy spectrum with an accuracy of about 5% up 
to a; g « 1.5. If we take this as a limit for the root mean 
square of the relative momentum (q^)^/^ = 0.53 k in the 
many-body system, we can apply our many-body model 
up to aiK ~ 3. 



C. Functional Variation and the Extremum 
Condition 

The ground state density of a system is found by min- 
imizing the energy functional 

E[ki,...,ke]= Jd^x£[Ki,...,KE]{x). (57) 

for given particle numbers iV^. This constraint is im- 
plemented with help of a set of Lagrange multipliers /i^ , 
which are the chemical potentials of the different compo- 
nents. The Legendre transformed functional 



F[ki, ...,ks] 

= E[ki, . . . ,ks] - ^^jiV^ 



|d^xf [.!,..., ..H](f)-E04(^) ^^^^ 



d^x J-[ki, . . . , Ks] (x), 



has to be minimized by functional variation. A necessary 
but not sufficient condition for a set of local Fermi mo- 
menta {ki(x), . . . , Ks(x)} to minimize the transformed 
energy functional F[ki, . . . , k^] is stationarity, i.e., that 
the first variation of F[ki, . . . , k^^] with respect to all 
K(^{x) vanishes 

—F\Ki,...,Kf]=Q foraU ^. (59) 

This extremum condition is fulfilled if the derivative of 
the integrand J^[ki, . . . , ks\{x) with respect to all local 
Fermi momenta vanishes at each point x 



d 



dK^{x) 



J-[ni, . . . , k-b\{x) = for all x, ^. (60) 



Inserting expression ( p^ for the energy density and eval- 
uating the derivative results in the extremum condition 



(61) 



^y^^[34(a:) + 5.^(f)4(f)] 



for all X and each ^. This is a coupled set of 
S polynomial equations for the local Fermi momenta 
{ki{x), . . . , ke:(x)} at some given point x. Note that the 
trivial solutions k^(x) = were separated already. Any 
real solution of the extremum condition ( |6l| ) corresponds 
to a stationary point of the energy functional. In general 
one has to check explicitly whether they correspond to 
a minimum of the energy functional or whether they are 
maxima or saddle points. 

All following investigations on the structure and sta- 
bility of degenerate Fermi gases and on the influence of 
s- and p-wave interactions are based on the extremum 
condition (^ij). Many physical conclusions can be drawn 
from its algebraic structure already. We will discuss these 
questions in detail for the one- and two-component Fermi 
gas in section ^ and ^, respectively. 



IV. SINGLE-COMPONENT FERMI GAS 

As a first application of the formalism developed in the 
preceding section we study the properties of a degenerate 
single-component Fermi gas. 



A. Effect of the p-Wave Interaction 

The energy density of the interacting multi-component 
Fermi gas (^) reduces for the single-component system 
to the form 



£[n]{x)^^U{x)K\x) 
67r^ 



1 



207r2 



m 



K^{x) 



(62) 



307r3f 



where k{x) is the local Fermi momentum. The first term 
is the contribution of the trapping potential U (x) , the 
second term is the kinetic energy, and the third term 
describes the contribution of the p-wave interaction with 
a p-wave scattering length ai. As mentioned earlier the 
s-wave part of the interaction does not contribute in a 
system of identical fermions due to the Pauli principle. 
Therefore the p-wave part is the leading interaction term 
and there is no reason to neglect it from the outset. 

A first hint on the effects of the p-wave interaction are 
given by the density distributions for different values of 
the p-wave scattering length. The density distribution 
is obtained by solution of the extremum condition (61), 
which takes the simple form 



m[fj. - U{x)] = fiinix)] 

1 o 



(63) 



with /i(k) = -^k 



157r 



This 5th order polynomial equation for the local Fermi 
momentum k{x) is solved numerically for each point x. 
The chemical potential /i is adjusted such that particle 
number (|5^) assumes the desired value. 
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FIG. 5: Density profile p(x) of a single-component Fermi 
gas of A*' = 10^ particles trapped in a spherical symmetric 
parabolic trap with oscillator length I. The solid curve shows 
the noninteracting gas ai/f = 0. Dotted curves correspond to 
repulsive p-wave interactions with a\/l = 0.03, 0.06, and 0.1 
(from top to bottom). The dashed curves show attractive p- 
wave interactions with ai/i = —0.03, —0.04 and —0.044 (top 
to bottom), respectively. 



Figure ra shows the resulting radial density profiles 
p{x) = K^x) / {'ott'^) for a single-component gas of = 
10^ particles in a spherical trap with oscillator length i 
for different p-wave scattering lengths ai. The oscillator 
length defines the fundamental length scale of the prob- 
lem and the parameter that determines the strength of 
the interaction is the ratio of the p-wave scattering length 
and oscillator length, ai/i. To increase the magnitude of 
this ratio experimentally one can either increase the mag- 
nitude of the scattering length or decrease the oscillator 
length. 

For a repulsive p-wave interaction, i.e., ai/^ > 0, of in- 
creasing strength (dotted curves) the density distribution 
flattens and expands radially compared to the noninter- 
acting system (solid line). For a ratio ai/i = 0.1 the cen- 
tral density has dropped to one half of the density of the 
noninteracting gas. With a typical experimental oscilla- 
tor length of ^ = 1/zm this ratio corresponds to a rather 
large scattering length of ai k, 2000aB, which neverthe- 
less may be within the range of experimental parameters 
[Q. For a tightly confining trap with I = O.l/xm a moder- 
ate scattering length of ai ~ 200aB is required to obtain 
the same ratio. 

For an attractive p-wave interaction, ai/£ < 0, the 
central density increases significantly with increasing in- 
teraction strength. If the central density exceeds a cer- 
tain value, i.e., if |ai/£| exceeds a critical value, then the 
extremum condition (p3) has no real solution any more. 
Physically this corresponds to a collapse of the dilute gas 
caused by the attractive mean-field that is generated by 
the p-wave interaction. We will discuss this question in 
detail in the following sections. 



The dependence of the density distribution on the 
p-wave scattering length as depicted in Figure ^ al- 
ready demonstrates that the p-wave interaction may have 
strong influence on the properties of degenerate Fermi 
gases. 



B. Mean-Field Instability: A Variational Picture 

To illustrate the origin and mechanism of the collapse 
of the metastable state of the trapped Fermi gas we utilize 
a simple variational picture. Assume a single-component 
Fermi gas of N particles in a spherical symmetric oscil- 
lator potential. The local Fermi momentum of the inter- 
acting system is parameterized by the analytic expression 
for the local Fermi momentum of the noninteracting sys- 
tem 



k(x) 



2{m) 



1/3 



for X < Xt 



(64) 



where the classical turning point X^ is treated as varia- 
tional parameter. By inserting this parameterization into 
the energy density (p3) and integrating we obtain a closed 
expression for the energy as function of the parameter Xt 



E{X, 



NX^ iV5/3 

— — r Ot 



with constant coefficients 



16m 



Ct = 



3(9/2) 
2m 



X^ 



1/3 



Ci- 



3 „3 



x^ 



(65) 



_ 8«(4/3)V3 

- 19257r2m ' ^'^^> 



Again the first term corresponds to the external poten- 
tials, the second to the kinetic energy, and the third term 
to the p-wave interaction. 

Figure ^ shows the dependence of the total energy ( |65| ) 
on the parameter Xt for a system of TV = 10^ particles 
with different p-wave scattering lengths. For attractive 
p-wave interactions, i.e., negative scattering length ai, 
the contribution of the interaction in (|6^) is negative. 



Due to its dependence this interaction contribution 
leads to a rapid drop of the energy for systems of de- 
creasing spatial extension and thus increasing density. 
At very high densities (small Xt) one formally ends up 
with states of negative energy, i.e., bound states. One 
should however keep in mind that the assumptions made 
for the construction of the Effective Contact Interaction 
are not valid in this high density regime any more. 

The not self-bound metastable state appears as local 
minimum at positive energies and low densities provided 
the p-wave attraction is sufficiently weak; the thick solid 
curve in Figure ^ shows an example. If the strength 
of the attractive p-wave interaction increases then the 
local minimum flattens and devolves to a saddle point 
(dashed curve) . From this particular interaction strength 
on the metastable low-density state does not exist any- 
more, only the true ground state of the system remains, 
which is usually a crystal. The system collapses if the 
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FIG. 6; Variational energy ( |65| ) of a trapped single-component 
Fermi gas with A'^ = 10^ as function of the parameter Xt. 
The curves show the noninteracting gas (thin solid), a\/l — 
-0.035 (solid), ai/£ = -0.051 (dashed), and a^/i = -0.065 
(dotted). The energies are given in units of the ground state 
energy Eo of the noninteracting gas. 




K [£-1] 

FIG. 7: Right hand side /i(k) of the extremum condition 
( |63[ ) as function of the Fermi momentum for a noninteracting 
singe-component gas (solid), with repulsive p-wave interaction 
ai/£ — 0.04 (dotted), and with attractive p-wave interaction 
ai/£ = —0.04 (dashed). The horizontal lines mark the respec- 
tive values of the chemical potentials for A'' = 10^ particles. 



barrier caused by the positive kinetic and the attractive 
mean-field energy vanishes. Since the mean-field attrac- 
tion grows with increasing density the system is unstable 
and collapses towards a high-density configuration. 



C. Mean-Field Instability: Stability Conditions 

Based on the extremum condition (|63| ) we derive a set 
of analytic stability conditions that relate the maximum 
density of a metastable system with the p-wave scattering 
length. Part of this was already discussed in 

The mean-field instability of the system occurs at val- 
ues of the chemical potential /i and the scattering length 
ai, where the extremum condition ( |63| ) does not have 
a real solution any more. This is shown in a pictorial 
way in Figure |^ where the the right hand side /i (k) of 
the extremum condition is plotted as function of k for 
different ai/£. The solution of the extremum condition 
at some specific point x is given by the value of k at 
which the respective curve reaches the value m[fi-'U{x)]. 
In the minimum of the trapping potential (we assume 
U{x) — in the minimum) the solution is given by the 
point where /i(k) reaches the value m/i. By moving 
towards the outer regions of the trap m[/Lt — U{x)] de- 
creases and one scans /i(k) down along ordinate until 
one reaches m[fi — U{x)] = 0, i.e., the classical turning 
point. 

For repulsive p-wave interactions (dotted curve) the 
r.h.s. of ( |63| ) is a monotonic growing function and so- 
lutions exist for arbitrary values of m[/i — U{x)]. If the 
scattering length ai is negative (dashed curve) then /i(k) 



chemical potential ^„ 



2ai 



3(37r 



,3/2 



(67) 



For values of the chemical potential /i > /imax no solution 
of the extremum condition exists, i.e., there is no meta- 
stable low-density state. Equivalently only solutions with 
local Fermi momenta below Kmax correspond to minima 
of the energy functional; those above Kmax (gray segment 
of the dashed curve) correspond to maxima of the energy. 
Thus we get a limiting condition for the local Fermi mo- 
mentum of the metastable state 



-aiK{x) < 



or in terms of the density 



-alp{x) < 



(68) 



(69) 



This is one form of the stability condition for the single- 
component Fermi gas. We note that this condition is 
completely independent of the trap geometry. As soon as 
the stability condition is violated somewhere in the trap, 
in general in the minimum of the trapping potential, the 
system will become unstable. 

For practical purposes we formulate a stability condi- 
tion in terms of the particle number N. The maximum 
particle number iVmax of the metastable degenerate Fermi 
gas is directly connected to the maximum chemical po- 
tential /imax- This relation is established numerically by 
solving the extremum condition for the maximum chemi- 
cal potential /imax and integrating over the resulting den- 
sity distribution to obtain the corresponding maximum 
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particle number (|5^). This is done for several scatte ring 
length ai assuming a deformed oscillator potential ( |49| ) 
with mean oscillator length £. Finally a parametrized 
form of the stability condition is fitted to this data. The 
parameterization is motivated by the nonintcracting gas, 
where the maximum local Fermi momentum is propor- 
tional to /£. Inserting this into the stability condi- 



tion (68) leads to the form 



C 



< 1 



c 



-2.246. 



(70) 



The parameter C is fitted to the numerical results, which 
are reproduced with a deviation far below 1%. Note that 
this condition is independent of the deformation of the 
harmonic oscillator trap ||l^ . 

For an interaction strength of ai/£ = —0.01 which 
corresponds to a scattering length oi « — 200aB for a 
trap with i = l/im the maximum particle number is 
-^max — 7.8 X 10^. This seems to be out of the range 
of present experiments. Nevertheless if we increase the 
strength of the p-wave attraction to ai/£ = —0.1 then 
the maximum particle number drops to A^max = 7800. 
Experimentally this could be achieved by utilizing a p- 
wave Feshbach resonance to increase the p-wave scatter- 
ing length to oi ~ — 2000aB as proposed by J. Bohn |^ 
for the ''"K system. 



V. TWO-COMPONENT FERMI GAS 

As second application we consider the degenerate two- 
component Fermi gas. 

A. Interplay between s- and p-Wave Interaction 

The general energy density of a trapped multi- 
component Fermi gas in Thomas-Fermi approximation 
( |5^) takes for the two-component system the following 
form 



£[ki,K2]{x) 



1 

6^ 
1 

207r2m 



[Ui{x)kI{x) + U2{x)4{x)] 

[4{x) + 4ix)] 



9TT^m 

-^[.?(f) + .l(f)] 
307r-^TO 

aj + bp 



(71) 



[kI{x)4{x) + 4{x)4ix)], 



where ki (x) and K2 (x) denote the local Fermi momenta of 
the two components. In contrast to the single-component 
system, both, s- and p-wave terms of the Effective Con- 
tact Interaction contribute. The s-wave interaction acts 
only between particles of different species and generates a 
contribution proportional to the product of the densities 



of both components in the energy density ([71|) . The p- 
wave term acts between particles of different components 
as well as between particles of the same species. For rea- 
sons of simplicity we assume the same p-wave scattering 
length oi for these different interactions. 

Including the constraint of given particle numbers A^i 
and of the two components with help of the chemi- 
cal potentials fii and fi2 (see section III C ) leads to the 
transformed energy density 



T[ki,K2]{x) = £[Ki,K2]ix) 



^1 ..3 



67r^ 



kHx) 



67r2 



4{x). 



(72) 



Functional variation of the transformed energy functional 
leads to the extremum condition. For the two-component 
system the general form ( |6l| ) reduces to a coupled set of 
two polynomial equations 



m[^i — Ui{x)] — -kI{x) + -^kI(^) + Tf'"'^'i(^) 



3n 



■ bo 



307r 



3k2{x) + 5k1{x)K2{x)\ , 



(73) 



where the second equation is generated by the exchange 
Ki{x) K2{x) and [/ii — Ui{x)] [fi2 — U2{x)]. Trivial 
solutions with ki{x) = and K2(a;) = 0, resp., are already 
separated in this expression. 

These coupled equations have a great variety of solu- 
tions. In order to show the generic phenomena of the 
two-component system without too many parameters we 
restrict ourselves to equal numbers of particles in both 
components iV = TVi = iV2 as well as trapping po- 
tentials that differ only by an additive constant, thus 
fi-U{x) ^ fii- Ui{x) = Ai2 - U2{x). 

We will concentrate the further studies on the stabil- 
ity of the degenerate two-component Fermi gas against 
mean-field collapse. For this phenomenon solutions with 
identical local Fermi momenta for both components, 
k(x) = Ki{x) = K2(x), are relevant. Under this assump- 
tion the extremum condition reduces to a single equation 



-U{x)] = /2[k(x)] 

with f2{K)^-K^ + -^K^ + -^H^. 

2 6Tr OTT 



(74) 



For simplicity we introduce a modified p-wave scattering 
length [21] 



ai = ai + bo/3, 



(75) 



which contains the s-wave effective volume parameter. 
In the following we will discuss the properties of the two- 
component Fermi gas as function of the s-wave and the 
modified p-wave scattering length. 

For other phenomena, like the separation of the 
two components due to repulsive interactions, different 
classes of solutions become important. We will discuss 
these in a forthcoming publication. 
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B. Mean-Field Instability: Stability Conditions 

The stability of the two-component Fermi gas under 
the influence of s- and p-wave interactions can be inves- 
tigated with tools similar to the single-component case. 
Here the solutions with identical Fermi momenta for both 
components are of interest. 

Similar to the single-component case the right hand 
side /2(k) of the extremum condition (74) may exhibit a 
maximum if the s-wave or the p-wave scattering length 
is negative. Thus the density and particle number of 
the metastable low-density state may be limited. For a 
detailed analysis one has to look at all possible combi- 
nations of signs of the s- and p-wave scattering lengths 
separately: 

ao > 0,ai > : For a purely repulsive interaction /2(k) 
is a monotonic growing function and no mean-field 
induced collapse occurs. 

ao < 0,ai < : For purely attractive interactions /2(k) 
shows a maximum; thus the density of the meta- 
stable low-density state is limited. 

ao > 0, &! < : The negative contribution of the p-wave 
interaction dominates /2(k) at high densities and 
generates a maximum, i.e., the mean-field induced 
collapse can occur even if the s-wave interaction is 
repulsive. 

ao < 0,ai > : It depends on the relative strength of 
the s- and p-wave interaction whether the r.h.s. of 
the extremum condition has a local maximum or 
grows monotonically. 

Especially the stabihty in the last case depends on a sub- 
tle competition between s- and p-wave interaction. More- 
over it shows some completely new phenomena, which 
will be discussed in the following section. 

For those cases where /2(k) has a maximum the value 
of the local Fermi momentum Kmax at the maximum is 
given by the equation 



-aoKn 



2[aiKn 



2' 



(76) 



Again Kmax is an upper limit for the local Fermi mo- 
menta, which can occur for a metastable low-density 
state of the two-component gas. Thus we can formulate 
the stability condition 

-aaK{x) - 2[ai«;(f)]3 < | (77) 

or equivalently in terms of the density 

-[67r2 aijp(x)]i/3 ~ 127r2 dlp{x) < |. (78) 

If these stability conditions are violated than no metasta- 
ble low-density state exists for the two-component Fermi 
gas. For a pure s-wave interaction (ai = 0) the stability 
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FIG. 8: Contour plot of the logarithm of the maximum par- 
ticle number, logigA'max, as function of the s- and p-wave 
scattering lengths for a two-component Fermi gas in a har- 
monic oscillator potential with mean oscillator length £. Se- 
lected contours are labeled with the corresponding value of 
logj^pA'max- In the white area at positive p-wave scattering 
length no collapse can occur, i.e., the maximum particle num- 
ber is infinity. 



TABLE I: Parameters of the fitted stability condition for 
the two-component Fermi gas for different interaction types. 



interaction type 


Co 


Ci 


Coi 


n 


ao < 0, (Ji < 


-1.835 


-2.570 


0.656 


1 


ao > 0, ai < 


-1.378 


-2.570 


1.360 


1 


ao < 0, fii > 


-1.835 


-1.940 


2.246 


3 



condition ( |77| ) reduces to the form —aoK{x) < 7r/2, which 
was obtained earlier by M. Houbiers et. al. Com- 
pared to this simple form the inclusion of the p-wave 
interaction reveals several new effects. 

Before we discuss the structure of ( |77| ) we formulate 
an equivalent stability condition in terms of the number 
of particles N — Ni — N2 oi each component. For given 
values of the two scattering lengths the maximum local 
Fermi momentum and the maximum chemical potential 
is calculated. From the solution of the extremum condi- 
tion (|7^) for these parameters the corresponding maxi- 
mum particle number is determined. This numerical data 
is fitted by a suitable parameterization of the stability 
condition in terms of the particle number and the two 
scattering lengths ao/i and ai/£ 



(79) 



Cn+l 
01 



< 1. 



This parameterization is constructed in analogy to the 
single-component case (|7^); the additional cross-term is 
necessary to achieve a similar accuracy with typical devi- 
ations below 1%. The parameters Cq, Ci, and Coi have 
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to be fitted for each combination of signs of the two scat- 
tering length separately. The value of n is not included in 
the fitting procedure but chosen by hand. The resulting 
values are summarized in Table |. 

Figure |^ illustrates the dependence of the maximum 
particle number resulting from this stability condition 
on the s- and p-wave scattering lengths. The contour 
plot shows the logarithm of the maximum particle num- 
ber for each component as function of ao/£ and di/£. 
The first gross observation is that attractive s- and p- 
wave interactions with similar scattering lengths set sim- 
ilar restrictions to the stability of the two-component 
Fermi gas. For example, a pure s-wave interaction with 
ao/£ ~ —0.05 leads to a maximum particle number of 
Nmax ~ 1.7 X 10^. In comparison a pure p-wave interac- 
tion with same scattering length ai/£ = —0.05 causes 
a collapse at even lower particle numbers of A^max ~ 
2.2 X 10^ 

If both interaction parts are attractive they cooper- 
ate and cause an instability at lower particle numbers or 
densities. 

If the interaction is attractive in one partial wave and 
repulsive in the other than the repulsive part leads to 
a stabilization, i.e., it increases the maximum particle 
number. Here a significant difference between s- and p- 
wave interactions arises: The stabilization caused by a 
repulsive s-wave interaction is rather weak. Compared to 
a pure p-wave interaction with ai/£ = —0.05 the presence 
of a s-wave repulsion of the same magnitude ao/£ = 0.05 
increases the maximum particle number only from 2.2 x 
10^ to 8.9 X 10^. In the opposite case of an attractive 
s-wave interaction a p-wave repulsion of same magnitude 
will always lead to an absolute stabilization, i.e., there is 
no collapse for arbitrary large particle number despite of 
the s-wave attraction. We will study these special effects 
in detail in the following section. 

These results clearly demonstrate that it is necessary 
to include the p-wave interaction if the scattering length 
ai is roughly in the same order of magnitude as the s- 
wave scattering length. Even if the ratio of the scattering 
lengths, ai/flQ, is approximately 0.3 dramatic effects like 
the p-wave stabilization, which is discussed in the next 
section, can occur. As can be seen from Figure || the p- 
wave interaction may be neglected only if the ratio di / uq 
is smaller than 0.1. 



C. Mean-Field Instability: p-Wave Stabilization 

Several new phenomena occur due to the competition 
between an attractive s-wave (oq < 0) and a repulsive 
p-wave interaction (ai > 0). To understand the origin 
of these phenomena, which are a unique property of this 
type of interactions, we investigate the right hand side 
/2(k) of the extremum condition ([m]). Figure ^ depicts 
the dependence of /2(k) on the local Fermi momentum 
for a s-wave scattering length aa/£ = —0.05 and three 
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FIG. 9: Right hand side /2(k) of the extremum condition 
(^) as function of the Fermi momentum for attractive s-wave 
interaction with ao/£ = —0.05 and a repulsive p-wave inter- 
action with ai/£ — 0.014 (dotted curve), 0.015 (solid), and 
0.016 (dashed). The gray segments of the curves correspond 
to maxima of the energy density. 



slightly different positive p-wave scattering lengths in the 
range oq/^ = 0.014 ... 0.016. 

Due to the dominant K^-dependence any repulsive p- 
wave interaction causes /2 (k) to grow fast for large Fermi 
momenta. Thus the maximum is only local and does not 
determine necessarily the maximum Fermi momentum or 
chemical potential as in the cases with attractive or van- 
ishing p-wave interaction. If the p-wave repulsion is suf- 
ficiently strong the local maximum vanishes completely 
and /2(k) is a monotonically growing function. Then a 
solution of the extremum condition exists for any den- 
sity, chemical potential or particle number. An example 
is shown by the dashed curve in Figure ^. It can be seen 
from equation (|7j) that the local maximum disappears if 
the ratio of the two scattering lengths fulfills the condi- 
tion 



|ao| 



> 



3^2/3 



0.311. 



(80) 



If this condition is fulfilled the p-wave repulsion causes 
an absolute stabilization of the system against s-wave in- 
duced collapse. In this case the total mean-field contri- 
bution of the interactions is always repulsive and grows 
monotonically with density. 

Notice that the p-wave scattering length necessary for 
this stabilization is only approx. 1/3 of the modulus of 
the s-wave scattering length. Obviously the p-wave in- 
teraction may have drastic infiuence on the stability even 
if it is significantly weaker than the s-wave interaction in 
terms of scattering lengths. 

For weaker p-wave repulsions /2('*) still shows a local 
maximum and in addition a (local) minimum at larger 
Fermi momenta. Exanrples are shown by the solid and 
dotted curves in Figure 0. In this case the extremum con- 
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FIG. 10: Upper plot: /2(fi;) of the extremum condition ( |74[ ) 
as function of the Fermi momentum for an interaction with 
ao/£ = —0.05 and a\ll = 0.015. Lower plot: Distribution 
of local Fermi momenta k(x) for a spherical trap of oscillator 
length I for the three different chemical potentials marked in 
the upper panel. Solid curves show the equilibrium profiles, 
dotted curves show metastable configurations. 



dition has two separate branches that correspond to lo- 
cal minima of the energy density, which are separated by 
branch of local maxima (gray segments). The branch at 
lower Fermi momenta corresponds to the usual family of 
low-density solutions that were obtained with other types 
of interactions too. It ends up at the local maximum 
with Kmax givcu by @ and m^^ax = /2(Kmax)- The 
solution branch at higher Fermi momenta gives rise to 
a new family of high-density solutions, which are unique 
for this type of interaction. It is bounded from below by 
the local minimum (Kmin, Mmin) and raises up to arbitrary 
Fermi momenta and chemical potentials. 

For values of [/i — C/ {x)\ between /imin and /imax the ex- 
tremum condition has two solutions kiow and Khigh with 
/2(kiow) = /2 ('thigh), See Figure |l^. In equilibrium the 
one with lower free energy density (|7^) is realized. We 
define a chemical potential /^trans at which the free energy 



densities of both branches are equal; the value of /Xtrans 
can be determined numerically. Since we expect the so- 
lution k{x) to correspond to a minimum of the energy 
functional at each point for \\x — \J{x)\ < /itrans the 
low-density branch gives the equilibrium solution and for 
[/i — U{x)] > /itrans the high-dcnsity branch does. This 
gives rise to a Maxwell construction for the r.h.s. of the 
extremum condition (74) as illustrated in the upper plot 
of Figure M. The dotted parts of the lower and up- 



per branch in Figure 10 correspond to local minima and 
may occur as metastable states that eventually undergo a 
transition to the energetically lower equilibrium solution. 

The structure of the density distribution depends cru- 
cially on the value of the chemical potential m, i.e., the 
particle number. The upper plot of Figure |lO| shows the 
r.h.s. of the extremum condition ( |74| ) for an interac- 
tion with ao/e = -0.05 and ai/i = 0.015. The dashed 
horizontal lines mark three different chemical potentials. 
The lower plot shows the radial dependencies of the lo- 
cal Fermi momenta k{x) = (67r^p(x))^/^ for these three 
chemical potentials assuming a spherical trap with oscil- 
lator length £. For chemical potentials fi < /itrans — case 

(a) — the equilibrium solution is completely on the low- 
density branch and we obtain the usual smooth density 
profile. If /i > /itrans — casc (b) and (c) — then the equi- 
librium solution in the center of the trap is given by the 
high-density branch, while the solution for the outer re- 
gions of the trap is given by the low-density branch. Thus 
the equilibrium density profile shows a jump in density 
by typically one order of magnitude as one approaches 
the center of the trap. The location of the discontinuity 
is always given by the equation [^—U{x)] = /ttrans which 
reflects the mechanical equilibrium between the low and 
the high-density phase, i.e., the equality of the pressures 
p = —T[k]{x) at the boundary between the two phases. 

If the chemical potential is still below /imax — case 

(b) — then a solution with a smooth low-density profile 
all over the trap may exist as metastable state. Due to 
density fluctuations this state may undergo a transition 
to the energetically preferred equilibrium state, which 
includes the high-density phase. 

The physical origin of the high-density phase is quite 
intuitive: The usual stability is determined by a compe- 
tition between kinetic energy, which favors low densities, 
and attractive mean- field, which prefers higher densities. 
If the attractive mean-field becomes too strong, than the 
kinetic energy is not able to stabilize the system any- 
more and the mean-field collapse occurs. In case of the 
high-density phase the attractive mean-field generated 
by the s-wave interaction has already overcome the sta- 
bilizing effect of the kinetic energy. Nevertheless the col- 
lapse is prevented by the repulsive p-wave contribution 
which grows for higher densities faster than the s-wave 
attraction and inhibits a further increase of density. We 
call this new phenomenon p-wave stabilized high-density 
phase — in contrast to the low-density phase stabilized 
by the kinetic energy, which is still present in the periph- 
eral regions of the trap. 
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The situation discussed so far assumes a repulsive p- 
wave interaction that is sHghtly to weak to cause the 
absolute stabilization according to (|8^). If the p-wave 
strength is decreased further, then the values of /imin, 
Mmax, and /itrans also dccreasc. If the ratio of the p-wave 
and s-wave scattering lengths drops below the limit 



|aol 



< 



160 
7297r2 



0.281 



(81) 



then the chemical potential of the minimum is negative, 
Mmin < 0. An example is shown by the dotted curve in 
Figure I For even weaker p-wave interactions with 



|ao| 



< 0.274 



(82) 



fi can be negative. That means that the high-density so- 
lution forms a self-bound state independent of the trap- 
ping potential. Therefore as soon as the maximum chem- 
ical potential /Zmax is exceeded the gas collapses into a 
self-bound high-density state which is independent of the 
trap. 

We summarize the variety of structures that appear for 
interactions with attractive s-wave (oo < 0) and repulsive 
p-wave part (oi > 0) in the following list 

0.311 < ai/|ao|: the p-wave repulsion stabilizes the sys- 
tem for arbitrary densities and particle numbers 
with a smooth low-density profile. For ratios of 
the scattering lengths near the limit a smooth but 
significant increase of the central density occurs. 

0.274 < ai/|ao| < 0.311: For chemical potentials below 
Mmax or particle numbers below the corresponding 
maximum particle number ( [79| ) the usual low den- 
sity solution exists. Above these values the p-wave 
stabilized high-density phase appears, i.e., the den- 
sity in the central region of the trap is increased by 
typically one order of magnitude compared to the 
low-density profile in the outer regions. 

ai/|ao| < 0.274: Below /^,nax or Ny^nx a- stable solution 
with the regular low-density profile exists. Above 
the system collapses; the high-density phase is not 
stable anymore. 

This subtle dependence on the ratio of the scattering 
lengths is illustrated in Figure ^. The white region for 
"strong" repulsive p-wave interactions shows the domain 
of absolute stabilization. The solid line corresponds to 
the condition dSG) for absolute stabilization, the dashed 



line to condition (82) for the stability of the high-density 
phase. The small area between those lines represents 
the parameter region where the p-wave stabilized high- 
density phase exists if the maximum particle number is 
exceeded. 

To conclude this section we perform a Gedankenex- 
periment: Assume a two-component Fermi gas of A'^ = 
Ni = N2 = 60000 particles in each component trapped 




FIG. 11: Evolution of the density profile of a two-component 
system of = A''i = N2 = 60000 particles with p-wave scat- 
tering length ai/£ = 0.03 according to the Gedankenexperi- 
ment described in the text. The s-wave scattering length is 
tuned in the range ao/i = —0.095 (solid), —0.1 (dash-dotted), 
-0.101 (dashed), and -0.102 (dotted). 



in a spherical oscillator potential with £ — l^m. May 
the interaction be composed of a repulsive p-wave part 
with ai/i — 0.03 and an attractive s-wave compo- 
nent which can be tuned within a small range ap/-^ = 
—0.095, . . . , —0.102, e.g., by using a Feshbach resonance. 
Figure |ri| shows the evolution of the density profile of the 
Fermi gas if the strength of the attractive s-wave inter- 
action is increased slowly such that density fluctuations 
are negligible. For a^/i — —0.095 (solid curve) and —0.1 
(dash-dotted) we observe a smooth low density profile, 
where the central density increases slightly with increas- 
ing s-wave attraction. A dramatic change happens if the 
attraction is increased to a^/t = —0.101 (dashed). For 
this interaction strength the particle number of the sys- 
tem is already above the maximum particle number given 
by (^ and the high-density phase appears and occupies 
a rather large volume. Further increase of the s-wave 
attraction (dotted) causes a growth of the high-density 
phase. If the limit ao/€ = —0.11 is reached than part of 
the high-density component is self-bound and the system 
is expected to collapse. 



VI. SUMMARY AND CONCLUSIONS 

We formulated a simple and transparent model to de- 
scribe the structure and stability of degenerate multi- 
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component Fermi gases trapped in an external potential. 
In a first step we derived an Effective Contact Interaction 
(ECI) for all partial waves that reproduces the exact two- 
body energy spectrum when used in a mean-field model 
space. Including the s- and p-wave parts of the ECI 
we constructed the energy density of the inhomogenous 
Fermi gas in a mean-field calculation using the Thomas- 
Fermi approximation. By functional minimization of the 
energy we obtained a set of coupled polynomial equa- 
tions for the ground state density profile of the system. 
We showed that the combination of s- and p-wave inter- 
actions leads to a rich variety of phenomena in trapped 
degenerate Fermi gases. 

In the single-component system the p-wave part is the 
leading interaction term since s-wave scatterings are pro- 
hibited by the Pauli principle. Attractive p-wave interac- 
tions cause a mean-field instability of the one-component 
gas if a certain maximum density is exceeded. We de- 
rived explicit stability conditions in terms of the density 
or particle number and the p-wave scattering length. 

The interplay between s- and p-wave interactions leads 
to several new effects in the two-component Fermi gas. 
We discussed the dependence of the mean-field instability 
on the s- and p-wave scattering lengths and derived also 
for this case closed stability conditions. It turns out that 
attractive s-wave as well as attractive p-wave interactions 
can cause a mean-field collapse. In addition a repulsive 
interaction part leads to stabilization, i.e., an increase of 
the maximum possible density of the Fermi gas. 

Interactions with attractive s-wave and repulsive p- 
wave parts show several special properties. If the p-wave 
scattering length exceeds about 1/3 of the modulus of 
the s-wave scattering length the s-wave attraction is fully 
compensated and no mean-field collapse occurs anymore 
at high densities. In the transition region towards this 
absolute stabilization a distinct high-density phase may 
appear in the center of the trap, which is stabilized by 
the p-wave repulsion alone. 

We conclude that the p-wave interaction may have an 
important influence on the structure and stability of di- 
lute degenerate Fermi gases. Considering the simulta- 
neous s- and p-wave Feshbach resonances predicted for 
the system it can be foreseen that large values 
of p-wave scattering length will be available experimen- 



of absolute stabilization, which we discussed, to allow 
higher densities of the normal Fermi gas and thus in- 
crease the transition temperature to a superfluid state 
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tems with strong attractive s-wave interactions are fa- 
vor ed § ^. Here the mean-field instability limits the 
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sen repulsive p-wave interaction one could use the effect 
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